; mv idl.ps ~/Dropbox/Tex/paper.nssl/fig/theprof.eps
if !d.name eq 'PS' then begin
    device,xsize=22,ysize=12,yoffset=3
    !p.thick=2 & !x.thick=1 & !y.thick=1
end
!x.margin=[4,3] & !y.margin=[4,3]
!p.multi=[0,2,1]

dir1='/home/guerrero/eumag/nt-28-rf-movt/'
dir2='/home/guerrero/eumag/tac-28-rf-movt/'

file1a=dir1+'velocity.sav'
file1b=dir1+'khelicity.sav'
file1c=dir1+'mhelicity.sav'
restore,file1a
restore,file1b
restore,file1c
print,'itac = 0'
itac=0
@rhoprof
print,rh
m=(size(rurms2))[2]
rgood=where(r ge 0.72)
Hrho=-(1./(deriv(r,rh)/rh))*6.96e8
Hrhom=mean(Hrho[rgood])
print,'Hrho = ',Hrhom, ' m'
urms_av=sqrt(mean(rurms2[rgood,*]))
print,'urms = ',urms_av, ' m/s'
tau=Hrhom/urms_av
alphak=(total(khel[0:31,*],1)/double(m/2))/rr
alpham=sqrt(total(mhel[*,*]^2,1)/double(m))/rr
alpha=-tau*(alphak+alpham)/(3.)

plot,r[0:45],smooth(alpha[0:45],2),xrange=[0.61,0.96],charsize=1.5,xstyle=1, $
	yr=[-20,30],ystyle=1,thick=5,xtitle='!8r/R!6'+sunsymbol(),$
	title='c) !7a!6 [m/s]',/nodata,xticks=4
oplot,r[1:45],smooth(alpha[1:45],1),li=0,thick=10,color=cgColor('red')

file2a=dir2+'velocity.sav'
file2b=dir2+'khelicity.sav'
file2c=dir2+'mhelicity.sav'
restore,file2a
restore,file2b
restore,file2c
print,'itac = 1'
itac=1
@rhoprof
m=(size(rurms2))[2]
rgood=where(r ge 0.72)
Hrho=-(1./(deriv(r,rh)/rh))*6.96e8
Hrhom=mean(Hrho[rgood])
print,'Hrho = ',Hrhom, ' m'
urms_av=sqrt(mean(rurms2[rgood,*]))
print,'urms = ',urms_av, ' m/s'
tau=Hrhom/urms_av
print,'tau = ',tau, ' seg ,', tau/(365.*24.*3600.), ' yr'
alphak=(tau/3.)*(total(khel[0:31,*],1)/double(m/2))/rr
alpham=(tau/3.)*(total(mhel[0:31,*],1)/double(m/2))/rr
alpha=-(alphak-alpham)

oplot,r[1:l-3],smooth(alpha[1:l-3],2),li=0,thick=10,color=cgColor('blue')
oplot,[0.6,0.98],[0.,0.],li=1,thick=2
oplot,[0.68,0.68],[-40,40],li=1,thick=2
oplot,[0.77,0.77],[-40,40],li=1,thick=2

;;;;;;;;;;;;;;
;;;;;;;;;;;;;;
;;;;;;;;;;;;;;

file1a=dir1+'velocity.sav'
file1b=dir1+'khelicity.sav'
file1c=dir1+'mhelicity.sav'
restore,file1a
restore,file1b
print,'itac = 0'
itac=0
@rhoprof
print,rh
m=(size(rurms2))[2]
rgood=where(r ge 0.72)
Hrho=-(1./(deriv(r,rh)/rh))*6.96e8
Hrhom=mean(Hrho[rgood])
print,'Hrho = ',Hrhom, ' m'
urms_av=sqrt(mean(rurms2[rgood,*]))
print,'urms = ',urms_av, ' m/s'
tau=Hrhom/urms_av
print,'tau = ',tau, ' seg ,', tau/(365.*24.*3600.), ' yr'
eta_t = tau*(total(rurms2,2)/double(m))/3.

plot,r[0:45],eta_t[0:45],xrange=[0.61,0.96],charsize=1.5,xstyle=1, $
	/ylog,yr=[1e6,6e9],ystyle=1,thick=5,xtitle='!8r/R!6'+sunsymbol(),$
	title='d) !4g!6!Dt!N!6 [m!U2!N/s!N]',/nodata,xticks=4
oplot,r[0:45],eta_t[0:45],li=0,thick=10,color=cgColor('red')

file2a=dir2+'velocity.sav'
file2b=dir2+'khelicity.sav'
file2c=dir2+'mhelicity.sav'
restore,file2a
m=(size(rurms2))[2]
;
print,'itac = 1'
itac=1
@rhoprof
rgood=where(r ge 0.72)
Hrho=-(1./(deriv(r,rh)/rh))*6.96e8
Hrhom=mean(Hrho[rgood])
print,'Hrho = ',Hrhom, ' m'
urms_av=sqrt(mean(rurms2[rgood,*]))
print,'urms = ',urms_av, ' m/s'
tau=Hrhom/urms_av
print,'tau = ',tau, ' seg ,', tau/(365.*24.*3600.), ' yr'
eta_t = tau*(total(rurms2,2)/double(m))/3.

oplot,r[1:62],eta_t[1:62],li=0,thick=10,color=cgColor('blue')
;
oplot,[0.8,0.85],[1.3e8,1.3e8],li=0,thick=10,color=cgColor('red')
xyouts,0.86,1e8,'CZ02',charsize=2
oplot,[0.8,0.85],[5e7,5e7],li=0,thick=10,color=cgColor('blue')
xyouts,0.86,4.0e7,'RC02',charsize=2

oplot,[0.68,0.68],[1e6,1e10],li=1,thick=2
oplot,[0.77,0.77],[1e6,1e10],li=1,thick=2

;;;;;;
;;;;;;
;;;;;;


;oplot,[0.85,0.9],[340,340],li=0,thick=10,color=cgColor('red')
;xyouts,0.91,300,'CZ02',charsize=2
;oplot,[0.85,0.9],[120,120],li=2,thick=10,color=cgColor('blue')
;xyouts,0.91,100,'RC02',charsize=2

end
